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recommended  track  spacing  for  an  Archimedes  spiral  search  in  a  datum  search  problem 
was  studied.  To  find  a  recommended  track  spacing,  three  analytical  approaches  were 
explored.  Each  of  three  analytical  approaches  has  its  own  strengths  and  weaknesses.  This 
analysis  combined  three  analytical  functions  into  a  single  parameterized  expression.  To 
find  the  best-fit  parameters  maximizing  probability  of  detection,  a  simulation  experiment 
with  a  NOLH  (Nearly  Orthogonal  Latin  Hypercube)  design  was  used. 
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I.  INTRODUCTION 


A.  BACKGROUND 

The  Republic  of  Korea  (ROK)  is  surrounded  by  several  powerful  countries  and 
also  by  water  on  three  sides.  For  several  decades,  the  ROK  has  faced  an  increasing  threat 
of  hostile  submarines.  In  response,  the  government  and  military  have  improved  their  anti¬ 
submarine  warfare  (ASW)  capabilities.  Weapon  systems  such  as  the  P-3C  maritime 
patrol  aircraft  (MPA)  have  been  deployed. 

These  ASW  systems  are  important  for  preventing  terrorist  invasions  by  submarine 
and  to  protect  the  sea  lanes.  This  thesis  will  develop  an  ASW  search  tactic  for  an  MPA 
aircraft  searching  for  a  patrolling  submarine. 

B.  SCENARIO 

The  periscope,  communications  antenna,  or  snorkel  of  a  target  submarine  is 
detected  by  the  long-range,  surface  search  ASW  radar  on  an  MPA  aircraft.  The  initial 
radar  detection  provides  an  estimated  target  position,  called  datum.  The  aircraft  then  flies 
directly  to  datum  and  conducts  a  spiral-out  ASW  search  for  the  target.  The  target 
submarine  submerges  below  periscope  depth  immediately  after  radar  detection  and 
continues  its  underwater  patrol. 


Figure  1.  ASW  Search  Scenario 
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C.  RESEARCH  QUESTION 

Referring  to  Figure  1,  the  problem  parameters  are: 

U  :  target  speed  ( kt ) 

V  :  aircraft  speed  ( kt  ) 

R  :  range  at  which  the  aircraft  can  detect  the  submarine  (mn) 
t  :  time  required  for  aircraft  to  fly  from  where  detection  occurred  to  datum 

(hr) 

A  :  the  Poisson  rate  at  which  the  target  changes  course  ( hr  1 ) 

The  main  research  question  investigated  here  is  as  follows:  Given  problem 
parameters  U,  V,  R,  r  ,  and  A  ,  what  is  the  aircraft  search  track  spacing  S*  which 
maximizes  the  probability  of  detection  of  the  submarine  by  the  ASW  aircrafts? 

D.  STRUCTURE  OF  THESIS 

Chapters  II  and  III  discuss  the  mathematical  model  of  a  target  submarine  and  a 
searcher  aircraft  which  conducts  an  ASW  search  operation.  Chapter  IV  develops  three 
analytical  approaches  to  suggest  good  track  spacings  for  the  spiral-out  search.  Chapter  V 
develops  the  design  of  the  simulation  experiment  and  the  analysis  for  constructing  a 
meta-model  to  find  the  recommended  track  spacing.  Chapter  VI  discusses  the  simulation 
results.  Chapter  VII  provides  conclusions  and  recommendations  for  further  research. 
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II.  SEARCHER  MOTION  MODEL 


A.  ARCHIMEDES  SPIRAL  SEARCH 


Initial  target  detection 


Aircraft  arrives  at  datum  and 
ASW  search  begins 


-r 


t  =  0 


Figure  2.  Search  Timeline 


The  target  is  initially  detected  at  time  t  =  —r.  At  time  0,  the  searching  aircraft 
arrives  at  datum  and  begins  a  spiral-out  search  with  constant  track  spacing  S  and  speed  V. 
The  aircraft  carries  a  detector  with  certain  detection  range  R. 

The  aircraft  track  will  follow  an  Archimedes  Spiral  [1],  For  such  a  search  path, 
the  radial  distance  from  datum  in  terms  of  the  angular  position  6  in  radians  is 


m= 


0S_ 

2n 


So  when  0  increases  by  2 n  (one  revolution),  the  radial  distance  increases  by  S. 


In  Chapter  IV,  it  is  shown  that  at  time  t  >  0 ,  the  approximate  radial  position  at  time  t  of 
a  searcher  with  speed  V  is 


r(t)  = 
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III.  TARGET  MOTION  MODEL 


A.  RANDOM  TOUR 

This  model  assumes  the  target’s  motion  is  a  random  tour  [2].  The  target’s  motion 
is  assumed  to  be  independent  of  the  search,  and  not  specifically  directed  toward  escape. 
The  target  may  move  a  considerable  distance  from  the  point  of  initial  detection  (datum) 
as  the  search  for  it  proceeds.  The  target  motion  can  be  beneficial  for  the  searcher  in  that  it 
increases  relative  speed,  but  it  can  also  be  harmful  as  it  increases  the  size  of  the  area  in 
which  the  target  can  be  located. 


Figure  3.  Description  of  Random  Tour  Target 

B.  PROBLEM  PARAMETERS 

Figure  3  shows  an  example  of  a  random  tour  target  motion  with  X  of  .5  course 
changes/hr,  and  U  of  20  nm/hr. 

1.  U  (Target  Speed) 

Target  speed  U  is  assumed  to  be  constant. 
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2. 


A,  (Course  Change  Rate) 


The  target  changes  its  direction  at  times  determined  by  a  Poisson  Process  [3]  with 
rate /l .  For  instance,  if  A  =2  /  hour,  the  target  changes  its  course  two  times  per  hour  on 
the  average,  and  the  time  between  course  changes  is  exponentially  distributed  with  mean 
1/2  hour. 


3.  <j  (Observation  Error) 

Immediately  after  the  initial  radar  detection  at  time  t  =  -r,  the  target’s  location 
distribution  is  circular  bivariate  normal  with  equal  means  (0,  0)  and  variances  (cr2,cr2  ) . 

In  this  thesis,  we  will  assume  o  =  1  nm,  consistent  with  a  long-range  radar  detection  of  a 
small  target. 

C.  NORMAL  APPROXIMATION  FOR  THE  PURE  RANDOM  TOUR 
TARGET 

For  large  enough  time  t ,  the  Central  Limit  Theorem  [4]  guarantees  that  the  target 
position  will  be  approximately  circular  bivariate  normal. 


XJ  :  target  speed 


T  :  duration  of  1st  leg 

(exponentially  distributed  with 

mean  y^) 

9  :  Target’s  course  during  1st  leg 
(uniformly  distributed  between 
0  -  2;r  radians) 

D  :  length  of  1st  leg 

(exponentially  distributed  with 

mean  ) 


Figure  4.  The  Initial  Random  Tour  Leg 
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To  see  this,  consider  the  following  development.  Since  target’s  course  change  rate  is  A , 
the  time  to  the  first  course  change,  T ,  is  an  exponentially  distributed  random  variable 
with  rate  A  .  We  write  this  as  T  ~  cxp(A) .  So, 

P(D<d)  =  P(UT<d ) 

=  p(t<%) 

=  l-exp 

Thus,  D-exp^-j^j.  That  is  D,  the  distance  travelled  by  the  target  during  the 
first  leg,  is  an  exponentially  distributed  random  variable  with  mean  . 

Now,  we  want  to  find  the  expected  value  and  variance  of  the  target  position  after 
one  leg.  Each  new  course  6  is  uniformly  distributed.  So, 

6  and 

X  =  DcosO. 

Now  by  symmetry  of  the  cosine  function, 

£(X)  =  0. 


So, 


V(X)  =  E(x'-)-(E(X)f 

=4V) 


Now  conditioning  on  the  target  course  6 , 

In 

e{x2)=  J  E[x2\d  =  T)fg{T)dT 

r=0 
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And  by  the  symmetry  of  cos2  (r) , 


|  cos2  zE^D2^dz . 

T= 0 


And  from  the  properties  of  the  exponential  distribution. 


cos2  zdz 


Now  by  x,y  symmetry, 


is  (7)  =  0 ,  and 


Let  XTot  be  the  x  coordinate  of  the  target  position  at  time  t . 

Then, 

XTol  =  Xl  +  X 2  ■+ - b  XN ,  where  N  ~  Poisson  (At) . 

Since  variances  add  for  independent  random  variables, 

V(xn,\ N  =  n)  =  fjV(X,) 


Now,  conditioning  on  N , 
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00 


^rj  =  l>(( U/A2)p(N  =  n) 


Thus, 

nxn,h{uyA)‘, 

v^,Huyy 

r(V,,)  =  0,:iiul 

E(YTc,)  =  0. 


And  by  the  Central  Limit  Theorem,  after  a  sufficiently  large  number  of  turns,  we  can 


approximate  the  target’s  position  with  a  circular,  bivariate  nonnal  random  variable  with 


mean  (0,0)  and  variances  in  the  x  and  y  directions  of 
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IV.  ANALYTICAL  APPRO ACHS  TO  FIND  GOOD  TRACK 

SPACINGS 


Let  A(t)  be  the  area  enclosed  by  the  track  by  time  t.  So, 

A(t)«SVt. 

Assuming  that  the  area  enclosed  is  approximately  circular, 

A{t^  ~  nr 2  (t)  . 

Equating  these  expressions  gives 

r(t)*JsVt/n,t>0.  (1) 

As  illustrated  in  Figure  5,  we  can  use  these  expressions  to  find  the  track  spacing 
( S* ),  such  that  r  (t)  is  just  tangent  to  U(t  +  r)-R  at  some  time  ( t  ). 

Referring  to  Figure  5,  at  time  t  =  t  , 

U  (t  +  r)-  R  =  r  (t)  =  yJs*Vt/n  . 


So, 


—  =  (u(t  +  v)-R)2 

=  {Ut  +  (Ur-R)f 
=  U2t 2  +  2 U  {Ur  -R)t  +  (Ur-R)2 . 
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Figure  5.  Calculation  of  S* 


Collecting  terms, 

UY  +  (lU (Ut - R) - SV/x)t  +  (Ut -R)2  =  0  . 

a  v - v - '  v 

b  c 

This  is  a  quadratic  equation  in  t .  There  can  be  0,  1  or  2  solutions.  We  are  looking 
for  the  S*  that  produces  one  solution.  From  the  quadratic  formula,  4 ac  =  b2  yields  one 
solution. 


So, 

4 U2  (Ut  -  R )2  =  (2 U  (i Ut  -  R)  -(SV) /xf 

.2  4 U(Ut-R)SV  ,  ,  x2 

=  4 U2  (Ut  -R) - k - l - +  (SV/7T )  . 

7T 

One  solution  is  S  =  0  .  Assuming  that  S  ^  0 ,  and  then  dividing  both  sides  by  S,  gives 

-4 ttU(Ut-R)V  +  SV2  _Q 
n  n1 


Solving  for  S  gives 
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(2) 


4n  U(Ut-R) 


S*  is  the  largest  track  spacing  which  keeps  all  search  effort  on  the  furthest-on- 


B.  SB 


At  time  t  =  —r  (initial  detection),  the  target’s  location  is  circular  normal  with 
means  (0,  0)  and  variances  (cr2,cr2)  .  The  target  then  begins  a  Random  Tour  with 
parameters  A  and  U  (target’s  speed). 

By  the  Central  Limit  Theorem,  the  target  position  at  time  t>—r  (for  sufficiently 
large  t)  is  still  circular  normal  with  mean  (0,  0),  but  with  variances  increased  from  cr2  (at 
,  U2(t  +  T ) 

time  t  =  -t  )  to  cr  -\ - - - -  (at  time  t>-r).  And  since  the  radial  position  of  a  circular 

A 

bivariate  normal  random  variable  with  variances  s 2  is  Rayleigh  distributed  with  mode  s, 
we  know  that  the  radial  position  at  time  t,  R(t) ,  is  approximately  Rayleigh  with  mode 


2  U2(t  +  Z) 
a+ - ^ 
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fir) 


distance  from  datum,  r 

Figure  6.  Density  Function  of  Radial  Position  of  Target  at  Time  t ,  Assuming  a  Normal 
Approximation  to  the  Random  Tour 

So,  a  spiral  search  which  maintains  a  radial  distance  at  search  time  t  of 

,  t  >  0  should  be  a  good  search.  The  justification  is  that  this  radial 

distance  is  the  mode  of  the  approximate  radial  distance  density  function. 

The  Archimedes  spiral  search  which  provides  the  tightest  support  below  for  this 
function  is 


as  illustrated  in  Figure  7. 
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Figure  7.  Calculation  of  SRay 


Setting 


U2t 

-  equal  to 

A 


(from  Equation  (1))  and  solving  for  S  yields 


SVt 

n 


U2t 

A 


so 


U27I 


VA  ' 


(3) 


SRay  should  be  a  good  track  spacing  for  a  random  tour  target  that  is  approximately 
normally  distributed.  However,  a  small  A  can  make  SRay  unreasonably  large  because  the 
nonnal  approximation  may  not  apply  in  that  case.  For  instance,  when 

1  0”  X  7t 

U  =  10kts,F  =  200kts,  A  =  0.0001/hour  ,  SRav  is  - «  15,700  miles,  which  is 

'  Ray  200x0.0001 

too  large. 
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c. 


S  =  2R 


When  the  target  is  slow  or  the  searcher  is  fast  (i.e.,  is  small),  then  an 

exhaustive  search  (track  spacing  of  2 R)  should  be  effective.  2 R  is  the  largest  track 
spacing  which  guarantees  no  gaps  in  search  coverage. 

D.  PROPOSED  EXPRESSION 

We  have  introduced  three  analytical  expressions  for  the  recommended  track 
spacing  :  SR  ,  S* ,  and  2 R.  Each  has  its  weaknesses  and  strengths.  In  particular, 

•  SRay  should  be  good  when  enough  time  has  elapsed  to  allow  the  Central 

Limit  Theorem  to  guarantee  that  the  target  distribution  is  roughly  bivariate 
nonnal.  But  as  we  have  seen,  when  X  is  small,  SRay  can  become 

arbitrarily  large.  Consequently,  it  should  be  bounded  above  by  another 
expression,  and  we  will  use  S*  for  this  purpose. 

•  S*  is  intended  to  keep  all  search  effort  within  the  furthest-on-disk.  As 
such,  S*  should  be  (approximately)  the  largest  track  spacing  considered. 
Otherwise  search  effort  would  be  placed  in  regions  where  the  target  can 
not  be  located. 

•  2 R  is  the  largest  track  spacing  which  guarantees  no  gaps  in  search 
coverage.  Thus  2 R  should  be  (approximately)  the  smallest  track  spacing 
considered  for  a  stationary  or  slow  target.  It  might,  however,  be  too  small 

when  R  is  small  or  the  speed  ratio  is  large. 

Considering  the  strengths  and  weaknesses  of  SR  ,  S* ,  and  2 R,  the  following 
expression  is  proposed  for  the  recommended  track  spacing: 

S  =  max{aiCmin{/?S^,S‘  +  y}}.  (4) 

Simulation  will  be  used  to  determine  best-fit  values  of  a ,  /? ,  and  y .  The  need 
for  y  was  discovered  when  conducting  the  simulation  experiments.  It  was  noted  that  as 
the  initial  observation  error  standard  deviation,  cr ,  increased,  the  best- fit  y  also  increased. 
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V.  DESIGN  OF  EXPERIMENT 


A.  NOLH  (NEARLY  ORTHOGONAL  LATIN  HYPERCUBE)  DESIGN 

Simulation  experiments  were  conducted  to  find  a  best-track  spacing  for  the 
Archimedes  spiral  search  pattern.  For  the  design  of  the  experiment,  a  full  factorial  design 
for  the  decision  variable  S  was  crossed  with  a  NOLH  [5]  design  for  the  other  five  factors. 
We  refer  to  the  factors  other  than  S  as  “Blocking  Factors”.  They  affect  the  outcomes,  and 
we  wish  to  understand  their  impact,  but  we  have  no  ability  to  change  them  in  a  real-world 
scenario  to  achieve  a  better  outcome.  To  construct  the  NOLH,  a  Microsoft  Excel 
spreadsheet  created  by  Professor  Susan  Sanchez  (Sanchez,  2005)  was  used.  The  number 
of  full  factorial  design  points  for  S  is  181  ranging  from  0  to  45nm  in  .25nm  increments, 
and  the  number  of  NOLH  design  points  is  65.  Appendix  A  is  an  extract  of  the  NOLH 
design  which  depicts  the  65  design  points.  The  65  design  points  are  close  to  orthogonal. 
(All  pairwise  correlation  values  for  design  columns  are  lower  than  0.05.)  The  two 
designs  are  crossed  to  create  the  following  design  table.  The  number  of  crossed  design 
points  is  65x181  =  11,765. 


DP 

U 

T 

V 

R 

A 

1 

14.656 

0.57 

153.91 

0.992 

0.25 

2 

19.109 

1.578 

116.41 

1.133 

0.688 

3 

17.922 

1.039 

242.97 

0.828 

0.594 

... 

65 

3.375 

1.016 

109.38 

0.711 

0.563 

NOLH  Design  for  Five  Blocking  Factors 


X 


DP 

5” 

1 

0 

2 

0.25 

... 

181 

45 

Full  Factorial 

Design  for  Decision  Variable  S 


DP 

U 

T 

V 

R 

k 

S=0 

S^O.25 

... 

S=45 

1 

14.656 

0.57 

153.91 

0.992 

0.25 

2 

19.109 

1.578 

116.41 

1.133 

0.688 

3 

17.922 

1.039 

242.97 

0.828 

0.594 

... 

65 

3.375 

1.016 

109.38 

0.711 

0.563 

Crossed  Design 


Figure  8.  Crossed  Design  for  Six  Factors 
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B.  THE  MAXIMUM  PROBABILITY  OF  DETECTION  (Pdbest)  FOR  EACH 
DESIGN  POINT 

Appendix  B  is  the  JAVA  code  for  the  search  simulation  model  used  for  data 
collection.  One  thousand  replications  are  used  for  each  of  the  11,765  crossed  design 
points.  By  running  the  simulation  model,  the  probability  of  detection  ( Pd  )  for  each 
design  point  was  estimated.  Table  1  and  Figure  9  show  some  of  the  Pd  results. 


Blocking 

Factor 

DP 

U 

T 

V 

R 

X 

5=0 

Pd 

5=0.25  ... 

5=45 

1 

14.656 

0.57 

153.91 

0.992 

0.25 

0 

0.054 

0.08 

2 

19.109 

1.578 

116.41 

1.133 

0.688 

0 

0.061 

0.078 

3 

17.922 

1.039 

242.97 

0.828 

0.594 

0 

0.089 

0.038 

64 

7.531 

1.813 

144.53 

1.156 

0.094 

0 

0.041 

0.069 

65 

3.375 

1.016 

109.38 

0.711 

0.563 

0 

0.397 

0.11 

Table  1 .  Probability  of  Detection  (  Pd  )  for  Each  Blocking  Factor  Design  Point 


Design  Point  1 


Figure  9.  Plot  of  Estimated  Simulation  Probabilities  of  Detection  for  Blocking  Factor 
Design  Point  1 
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Figure  9  shows  the  simulation  Pd  for  design  point  1  and  the  181  track  spacings 
(0-45nm).  We  can  now  determine  the  largest  Pd  ( Pdbest )  and  the  track  spacing  producing 
Pdbest  ( Shest ).  Example  results  are  given  in  Figure  9  and  Table  2. 


Blocking 

Factor 

DP 

U 

T 

V 

R 

/l 

^ best 

Pdbes, 

1 

14.656 

0.57 

153.91 

0.992 

0.25 

10.75 

0.758 

2 

19.109 

1.578 

116.41 

1.133 

0.688 

22.5 

0.116 

3 

17.922 

1.039 

242.97 

0.828 

0.594 

17.75 

0.39 

64 

7.531 

1.813 

144.53 

1.156 

0.094 

9.5 

0.835 

65 

3.375 

1.016 

109.38 

0.711 

0.563 

1.75 

0.943 

Table  2.  Table  for  Sbest  and  Pdbest  for  Each  Blocking  Factor  Design  Point 


C.  FINDING  THE  BEST  a,  (L  AND  y 

Using  the  simulation  data  of  Tables  1  and  2,  we  can  determine  the  values  of  a , 
(3 ,  and  y  in  Equation  (4)  which  produce  a  track  spacing  and  resulting  Pd  as  close  as 
possible  to  the  maximum  obtainable  Pd.  We  use  the  following  procedure: 

For  each  of  the  5F  possible  values  of  (<2,  /?,  /),  where  a, (3  and  y  vary  over 
{0,  .1,  .2,  5} 

For  each  of  the  65  blocking  factor  design  points  (dp) 

Determine  S (a,  f3,  y,  dp)  from  Equation  (4) 

Determine  Pd(s,  dp)  from  Table  1 


19 


Determine  Pdbest{ dp) 


from  Table  2 


*  best 


Determine  MOE,  (a,  p,  y,  dp)  =  ^yp^ 
Continue  for  each  of  the  design  points.  Then  compute 


MOE2 (a,  p,  y)  =  “tIXmOEj (a,  p,  y,  dp) 
vb5y  vdP 

Continue  for  next  (a,  p,  y)  to  find  by  enumeration 


(a*,  P* ,  y  j  =  argmax  MOE2(a,  p,  y) 

a,P,r 

For  each  (a,  p,  y)  triplet,  the  analysis  generates  65  values  of 
MOE,  (a,  p,  y,  dp)  -  one  for  each  of  the  blocking  factor  design  points.  Then 
MOE2(a,  P,  y)  is  the  average  of  these  65  values.  ( a  ,p* ,y are  the  values  which 
maximize  the  average  of  the  65  values  of  MOE, . 
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VI.  SIMULATION  RESULTS 


A.  MAIN  RESULT 

Using  the  enumeration  procedure  of  the  previous  section,  51'  (=  132,651)  values 
of  MOE2  («,  /?,  y)  were  computed.  The  largest  MOE2  was  .928,  and  the  maximizing 
arguments  were  (a*,/3*,y*^  =  (2,  2.7,  1.8). 


B.  VISUALIZATION  OF  MAIN  RESULT  USING  CONTOUR  PLOTS 


Contour  plots  were  constructed  to  help  visualize  how  MOE2  (a,  /?,  y)  varied 
over  the  range  of  its  arguments.  These  plots  reproduces  below  indicated  that  the  best 
values  of  (a,  j3,  y)  were  approximately (2,  2.7,  1.8) . 


Contour  Plot  for  Mean  Ratio 


Mean  Ratio 
■  <=  0.200 
■  <=  0.282 

■  <=  0.363 

■  <=  0.445 

■  <=  0.527 

■  <=  0.608 

■  <=  0.690 

■  <=  0.771 

■  <=  0.853 

■  <=  0.935 

■  >  1.000 


=  0.216 
=  0.298 
=  0.380 
=  0.461 
=  0.543 
=  0.624 
=  0.706 
=  0.788 
=  0.869 
=  0.951 


=  0.233 
=  0.314 
=  0.396 
=  0.478 
=  0.559 
=  0.641 
=  0.722 
=  0.804 
=  0.886 
=  0.967 


=  0.249 
=  0.331 
=  0.412 
=  0.494 
=  0.576 
=  0.657 
=  0.739 
=  0.820 
=  0.902 
=  0.984 


I  <=  0.265 
I  <=  0.347 
I  <=  0.429 
I  <=0.510 
I  <=  0.592 
I  <=  0.673 
I  <=  0.755 
I  <=  0.837 
I  <=  0.918 

I  <=  1.000 


Figure  10.  Contour  Plot  for  a  and  (5  given  y=  0 
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Figure  11. 


Figure  12. 


Contour  Plot  for  Mean  Ratio 


■  <=  0.200 

i  <=0.216 

<=  0.233 

<=  0.249 

<=  0.265 

■i  <=  0.282 

<=  0.298 

<=  0.314 

<=  0.331 

<=  0.347 

■  <=  0.363 

<=  0.380 

<=  0.396 

<=  0.412 

<=  0.429 

■  <=  0.445 

<=  0.461 

■1 

<=  0.478 

■1 

<=  0.494 

■■ 

<=  0.510 

■  <=  0.527 

<=  0.543 

<=  0.559 

<=  0.576 

■■ 

<=  0.592 

■  <=  0.608 

1  <=  0.624 

■■ 

<=  0.641 

■■ 

<=  0.657 

■■ 

<=  0.673 

■  <=  0.690 

!  <=  0.706 

<=  0.722 

<=  0.739 

<=  0.755 

N- 

N- 

O 

II 

i 

<=  0.788 

■1 

<=  0.804 

m 

<=  0.820 

<=  0.837 

■  <=  0.853 

SSI  <=  0.869 

<=  0.886 

<=  0.902 

<=  0.918 

■  <=  0.935 

■  >  1.000 

■  <=  0.951 

m 

<=  0.967 

<=  0.984 

<=1.000 

Contour  Plot  for  a  and  /?  given  y  =  1 


Mean  Ratio 
■  <=  0.200 
■  <=  0.282 

■  <=  0.363 
M  <=  0.445 
■i  <=  0.527 

■  <=  0.608 

■  <=  0.690 

■  <=  0.771 

■  <=  0.853 
■■  <=  0.935 

■  >  1.000 


=  0.216 
=  0.298 
=  0.380 
=  0.461 
=  0.543 
=  0.624 
=  0.706 
=  0.788 
=  0.869 
=  0.951 


■  <=  0.233 
<=  0.314 
<=  0.396 
<=  0.478 
<=  0.559 

■  <=  0.641 

■  <=  0.722 

■  <=  0.804 

■  <=  0.886 
■  <=  0.967 


=  0.249 
=  0.331 
=  0.412 
=  0.494 
=  0.576 
=  0.657 
=  0.739 
=  0.820 
=  0.902 
=  0.984 


■  <=  0.265 

■  <=  0.347 

■  <=  0.429 

■  <=  0.510 

■  <=  0.592 

■  <=  0.673 
I  <=  0.755 

■  <=  0.837 
I  <=  0.918 

I  <=1.000 


Contour  Plot  for  a  and  /?  given  y  =  2 


22 


Contour  Plot  for  Mean  Ratio 


Mean  Ratio 


■  <=  0.200 

■ 

<=  0.216 

■  <=  0.233 

■  <=  0.249 

■  <=  0.265 

■■  <=  0.282 

■ 

<=  0.298 

■  <=0.314 

■  <=  0.331 

■  <=  0.347 

■  <=  0.363 

H* 

<=  0.380 

<=  0.396 

■  <=  0.412 

■  <=  0.429 

■  <=  0.445 

<=  0.461 

<=  0.478 

■  <=  0.494 

■  <=0.510 

■  <=  0.527 

■ 

<=  0.543 

■  J  <=  0.559 

■  <=  0.576 

■  <=  0.592 

■  <=  0.608 

■ 

<=  0.624 

■  <=  0.641 

■  <=  0.657 

■  <=  0.673 

■  <=  0.690 

■ 

<=  0.706 

■  <=  0.722 

■  <=  0.739 

■  <=  0.755 

i-2 

o 

n 

i 

■ 

<=  0.788 

■  <=  0.804 

■  <=  0.820 

■  <=  0.837 

■  <=  0.853 

■i 

<=  0.869 

■  <=  0.886 

■  <=  0.902 

■  <=0.918 

■  <=  0.935 

<=  0.951 

■  <=  0.967 

■  <=  0.984 

■<=1.000 

>1.000 


Figure  13.  Contour  Plot  for  a  and  (5  given  y=  3 


Contour  Plot  for  Mean  Ratio 


Mean  Ratio 
■  <=  0.200 
■  <=  0.282 

■  <=  0.363 

■  <=  0.445 

■  <=  0.527 

■  <=  0.608 

■  <=  0.690 

■  <=  0.771 

■  <=  0.853 

■  <=  0.935 

■  >  1.000 


=  0.216 
=  0.298 
=  0.380 
=  0.461 
=  0.543 
=  0.624 
=  0.706 
=  0.788 
=  0.869 
=  0.951 


=  0.233 
=  0.314 
=  0.396 
=  0.478 
=  0.559 
=  0.641 
=  0.722 
=  0.804 
=  0.886 
=  0.967 


=  0.249 
=  0.331 
=  0.412 
=  0.494 
=  0.576 
=  0.657 
=  0.739 
=  0.820 
=  0.902 
=  0.984 


=  0.265 
=  0.347 
=  0.429 
=  0.510 
=  0.592 
=  0.673 
=  0.755 
=  0.837 
=  0.918 
=  1.000 


Figure  14.  Contour  Plot  for  a  and  (5  given  y=4 
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VII.  CONCLUSION  /  FUTURE  RESEARCH 


A.  CONCLUSION 

Since  the  military  threat  of  hostile  submarines  is  increasing,  the  need  for  effective 
ASW  operations  is  also  increasing.  In  this  thesis,  the  recommended  track  spacing  for  an 
Archimedes  spiral  search  in  a  datum  search  problem  was  studied.  This  thesis  has  sought  a 
meta-model  to  recommend  a  good  track  spacing.  This  analysis  combined  three  analytical 
functions  into  a  single  parameterized  expression.  To  find  the  best- fit  parameters 

maximizing  the  average  ratio  ^%),  ,  a  simulation  experiment  with  a  NOLH  (Nearly 

/  ‘  ^  best 

Orthogonal  Latin  Hypercube)  design  was  used.  The  complete  model  is 
S  =  max {2 R, min{2.7Srav , S*  + 1 .8} }  , 

where  R  is  detection  range,  S*  and  SRay  are  given  by  Equations  (2)  and  (3). 

B.  FUTURE  RESEARCH 

This  thesis  blended  different  techniques  involving  stochastic  modeling  and 
simulation,  design  of  experiments,  and  optimization  to  create  a  meta-model  for  track 
spacing  for  an  Archimedes  Spiral.  Suggested  future  research  includes: 

•  Good  search  patterns  with  multiple  targets  and  multiple  searchers; 

•  Selection  of  appropriate  spiral  search  pattern  in  datum  search  problem; 

•  Factor  analysis  for  probability  of  detection  in  ASW  datum  search; 

•  Recommended  track  spacing  when  the  other  parameters  are  unknown; 

•  Impact  on  observation  error  for  probability  of  detection. 
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APPENDIX  A.  NOLH  (NEARLY  ORTHOGONAL  LATIN 
HYPERCUBE)  DESIGN  FOR  FIVE  BLOCKING  FACTORS 


Low  level 

1 

0.5 

100 

0.5 

0 

High  level 

20 

2 

250 

2 

2 

decimals 

3 

3 

2 

3 

3 

factor  name 

U 

T 

V 

R 

A 

1 

14.656 

0.57 

153.91 

0.992 

0.25 

2 

19.109 

1.578 

116.41 

1.133 

0.688 

3 

17.922 

1.039 

242.97 

0.828 

0.594 

4 

13.172 

1.836 

207.81 

1.18 

0.125 

5 

18.516 

1.203 

128.13 

0.523 

0.188 

6 

11.094 

1.883 

135.16 

1.227 

0.313 

7 

15.547 

0.781 

179.69 

0.547 

0.5 

8 

16.438 

1.648 

238.28 

0.969 

0.75 

9 

14.063 

0.547 

102.34 

1.719 

0.813 

10 

19.406 

1.531 

172.66 

1.648 

0 

11 

10.797 

0.523 

245.31 

1.297 

0.406 

12 

19.703 

1.273 

203.13 

1.883 

0.344 

13 

11.391 

0.828 

149.22 

1.461 

0.781 

14 

15.844 

1.32 

163.28 

1.836 

0.531 

15 

11.984 

0.992 

217.19 

1.93 

0.875 

16 

13.766 

1.367 

182.03 

1.578 

0.281 

17 

17.328 

1.156 

158.59 

0.688 

1.344 

18 

11.688 

1.789 

139.84 

1.016 

1.156 

19 

14.953 

1.086 

184.38 

0.945 

1.563 

20 

12.578 

1.602 

228.91 

0.594 

1.094 

21 

16.734 

0.805 

160.94 

0.641 

1.969 

22 

16.141 

1.742 

100 

0.734 

1.031 

23 

20 

1.109 

219.53 

0.875 

1.938 

24 

12.281 

2 

224.22 

1.063 

1.375 

25 

12.875 

0.875 

114.06 

1.391 

1.625 

26 

15.25 

1.438 

123.44 

2 

1.531 

27 

17.031 

0.641 

212.5 

1.695 

1.063 

28 

18.219 

1.766 

193.75 

1.742 

1.844 

29 

18.813 

0.945 

118.75 

1.414 

1.281 

30 

14.359 

1.906 

151.56 

1.602 

1.781 

31 

13.469 

0.688 

205.47 

1.344 

1.906 

32 

17.625 

1.484 

240.63 

1.789 

1.438 

33 

10.5 

1.25 

175 

1.25 

1 

34 

6.344 

1.93 

196.09 

1.508 

1.75 

35 

1.891 

0.922 

233.59 

1.367 

1.313 

36 

3.078 

1.461 

107.03 

1.672 

1.406 

37 

7.828 

0.664 

142.19 

1.32 

1.875 

38 

2.484 

1.297 

221.88 

1.977 

1.813 

39 

9.906 

0.617 

214.84 

1.273 

1.688 

40 

5.453 

1.719 

170.31 

1.953 

1.5 

27 


Low  level 

1 

0.5 

100 

0.5 

0 

High  level 

20 

2 

250 

2 

2 

decimals 

3 

3 

2 

3 

3 

factor  name 

V 

T 

V 

R 

/l 

41 

4.563 

0.852 

111.72 

1.531 

1.25 

42 

6.938 

1.953 

247.66 

0.781 

1.188 

43 

1.594 

0.969 

177.34 

0.852 

2 

44 

10.203 

1.977 

104.69 

1.203 

1.594 

45 

1.297 

1.227 

146.88 

0.617 

1.656 

46 

9.609 

1.672 

200.78 

1.039 

1.219 

47 

5.156 

1.18 

186.72 

0.664 

1.469 

48 

9.016 

1.508 

132.81 

0.57 

1.125 

49 

7.234 

1.133 

167.97 

0.922 

1.719 

50 

3.672 

1.344 

191.41 

1.813 

0.656 

51 

9.313 

0.711 

210.16 

1.484 

0.844 

52 

6.047 

1.414 

165.63 

1.555 

0.438 

53 

8.422 

0.898 

121.09 

1.906 

0.906 

54 

4.266 

1.695 

189.06 

1.859 

0.031 

55 

4.859 

0.758 

250 

1.766 

0.969 

56 

1 

1.391 

130.47 

1.625 

0.063 

57 

8.719 

0.5 

125.78 

1.438 

0.625 

58 

8.125 

1.625 

235.94 

1.109 

0.375 

59 

5.75 

1.063 

226.56 

0.5 

0.469 

60 

3.969 

1.859 

137.5 

0.805 

0.938 

61 

2.781 

0.734 

156.25 

0.758 

0.156 

62 

2.188 

1.555 

231.25 

1.086 

0.719 

63 

6.641 

0.594 

198.44 

0.898 

0.219 

64 

7.531 

1.813 

144.53 

1.156 

0.094 

65 

3.375 

1.016 

109.38 

0.711 

0.563 
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APPENDIX  B.  JAVA  CODE  FOR  DATA  COLLECTION 


The  Java  code  within  this  appendix  collects  all  the  simulation  results  over  11765 

crossed  design  points. 

/* 

*  File:  Arch. java 

*  Created  on  Sep  15,  2007,  12:00  PM 
*/ 

j  kk 

*  Qauthor  Byungsoo  Son 

*  Class  to  collect  the  probability  of  detection 

*  for  the  11765  design  points 
*/ 

import  j ava . util . Random; 


public  class  Arch  { 


public  static  void  main (String [ ]  args)  { 

//set  input  values 

Random  generator  =  new  Random () ; 

//initializing  all  variables 

double  targetSpeed=0; 

double  tau=0; 

double  tMax  =5; 

double  sear cher Speed=0; 

double  detectionRange=0 ; 

double  trackSpacing=0 ; 

double  lambda=0; 

double  timeStep=0; 

double  xTarget=0; 

double  yTarget=0; 

double  sigma=0; 

double  xSearcher=0; 

double  ySearcher=0; 

double  distance=0; 

int  Detection=0; 

double  courseDirection; 

int  numberRelications  =  Integer .parse Int (args [ 0 ])  ; 
int  numberDetection=0  ; 

if  (args. length  ==  9) { 

sigma  =  Double .parseDouble (args [1 ])  ; 
targetSpeed  =  Double . parseDouble (args [2 ]) ; 
tau  =  Double. parseDouble (args [3] ) ; 
searcherSpeed  =  Double . parseDouble (args [ 4 ]) ; 
detectionRange  =  Double .parseDouble (args [ 5] ) ; 
lambda  =  Double .parseDouble (args [ 6] ) ; 
trackSpacing  =  Double .parseDouble (args [ 7 ]) ; 
timeStep  =  Double .parseDouble (args [ 8 ]) ; 

} 

else  { 

System. err. println ( "Supply  nine  arguments:  "  + 

"numberReplications (int) ,  sigma (double)  "  + 
"targetSpeed (double) ,  tau (double),  "  + 
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} 


"searcherSpeed (double) ,  "  + 

"detectionRange (double) ,  lambda (double) ,  "  + 
"trackSpacing (double)  "  + 

"and  timeStep (double) " ) ; 

System. exi t ( - 1 ) ; 

} 

//  replicate  n  times  given  the  other  parameters 
for(int  n=0 ; n<numberReplications ; n++) { 

xTarget  =  generator . nextGaussian () *sigma ; 
yTarget=  generator . nextGaussian ( ) *  sigma; 
xSearcher=0 ; 
ySearcher=0 ; 

Detection=0; 
cour seDirection  = 

generator . nextDouble ( ) *2  *Math . PI; 

//target  motion  model  for  time  -tau  to  0 
for (int  i=0;i< (tau/timeStep) ; i++) { 

if (generator . nextDouble ( ) <lambda*timeStep) { 
courseDirection  = 

generator . nextDouble ( ) *2*Math . PI; 

} 

xTarget  =  xTarget  +  Math. cos (courseDirection) 

* target Speed* timeStep; 

yTarget  =  yTarget  +  Math. sin (courseDirection) 

* tar get Speed* time Step; 

} 

//searcher  and  target  motion  model  for  time  0  to  tMax 
//If  the  distance  between  searcher  and  target  is  less 
//than  the  detection  range  of  sensor , detection  occurs 
for (int  j  =0 ; j < ( tMax/timeStep) ; j++) { 

if (generator . nextDouble  () <  lambda*timeStep)  { 

courseDirection  =  generator . nextDouble ( ) 
*2*Math . PI; 

} 

xTarget  =  xTarget  +  Math. cos (courseDirection) 

*  targetSpeed*timeStep; 
yTarget  =  yTarget  +  Math. sin (courseDirection) 

*  targetSpeed*timeStep; 

xSearcher  =  Math. sgrt (trackSpacing  *  searcherSpeed 

*  (timeStep  *  j) /Math. PI) 

*  Math . cos ( 2*Math . sgrt ( searcherSpeed 
*Math.  PI*  timeStep* j /trackSpacing) ) ; 

ySearcher  =  Math. sgrt (trackSpacing  *  searcherSpeed 

*  (timeStep  *  j) /Math. PI) 

*  Math . sin ( 2 *Math . sgrt ( searcherSpeed 
*Math. PI  *  timeStep* j /trackSpacing) ) 

distance  =  Math. sgrt (Math. pow(xTarget-xSearcher, 2) 

+  Math . pow( yTarget- ySearcher ,  2 ) )  ; 
if (distance  <  detectionRange) { 

Detection  =  Detection  +  1; 
break ; 

} 

} 

numberDetection  =  numberDetection  +  Detection; 

} 

System. out. println( (double) number Detect ion/ number Replications ) ; 
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